Undamped electrostatic plasma waves 
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Electrostatic waves in a collision-free unmagnetized plasma of electrons with fixed ions are investi- 
gated for electron equilibrium velocity distribution functions that deviate slightly from Maxwellian. 
Of interest are undamped waves that are the small amplitude limit of nonlinear excitations, such 
as electron acoustic waves (EAWs). A deviation consisting of a small plateau, a region with zero 
velocity derivative over a width that is a very small fraction of the electron thermal speed, is shown 
to give rise to new undamped modes, which here are named corner modes. The presence of the 
plateau turns off Landau damping and allows oscillations with phase speeds within the plateau. 
These undamped waves are obtained in a wide region of the (k, uj r ) plane (uj r being the real part of 
the wave frequency and k the wavenumber), away from the well-known 'thumb curve' for Langmuir 
waves and EAWs based on the Maxwellian. Results of nonlinear Vlasov-Poisson simulations that 
corroborate the existence of these modes are described. It is also shown that deviations caused by 
fattening the tail of the distribution shift roots off of the thumb curve toward lower fc-values and 
chopping the tail shifts them toward higher fc-values. In addition, a rule of thumb is obtained for 
assessing how the existence of a plateau shifts roots off of the thumb curve. Suggestions are made 
for interpreting experimental observations of electrostatic waves, such as recent ones in nonneutral 
plasmas. 



I. INTRODUCTION 

In his 1946 seminal paper [l| Landau demonstrated 
that electrostatic plasma waves of vanishing amplitude 
can be damped, due to their interaction with particles 
that stream with velocities close to the wave phase speed, 
Vf/,. For unmagnetized uniform plasmas, the wave damp- 
ing rate is generally proportional to the slope of the equi- 
librium distribution of particle velocities at v^. There- 
fore, for monotonically decreasing equilibrium velocity 
distribution functions (such as the usual Maxwellian) 
plasma waves are damped exponentially in time. 

Almost twenty years later O'Neil @ analyzed the ef- 
fects of nonlincarity on the propagation of plasma waves 
and found that the process of particle trapping in the 
wave potential well can inhibit Landau damping, by 
flattening the velocity distribution near the wave phase 
speed. 

In 1991 Holloway and Doming Q noted that certain 
nonlinear electrostatic oscillations can survive Landau 
damping even when their phase velocities are compa- 
rable to the electron thermal speed, Vth, due to the 
effect of particle trapping. They called these waves 
electron acoustic waves (EAWs), since in the range of 
small wavenumbers their dispersion relation is of the 
form to ~ 1.31kv t h- An EAW is in fact a Bernstein- 
Greene-Kruskal (BGK) mode [4j with a velocity distri- 
bution function that is effectively flat at the wave phase 
speed, due to trapping. In 1991, Demeio and Holloway [E[ 
performed Vlasov-Poisson simulations and provided evi- 
dence for the analytical results of Ref. [|[ . Also in Ref. 3 , 
the authors claimed that these undamped plasma oscilla- 
tions have no linear counterpart, since in their construc- 



tion the trapped particle distribution vanishes with wave 
amplitude and approaches a Maxwellian, for which the 
waves are heavily damped. 

The latter point was discussed in 1994 by Shadwick 
and Morrison [||, where it was pointed out that station- 
ary inflection point modes are the natural linear limit of 
the EAWs (BGK modes) of Ref. 0. A stationary inflec- 
tion point mode is an undamped mode with phase speed 
Vtj, at a point where the first two derivatives of a stable 
homogeneous equilibrium velocity distribution function 
vanish, and by the Penrose criterion this is necessary 
for the simultaneous existence of the mode and stabil- 
ity. The issue here is that the homogeneous equilibrium 
distribution function approached as the wave amplitude 
approaches zero is not unique, and there is no a priori 
reason it should be locally Maxwellian, since the struc- 
ture near is determined by the history of formation of 
the waves. Given a nonlinear BGK wave, there are many 
ways the limit of the vanishing of trapped particles can 
be taken and the result is clearly norm dependent. (See 
Ref. for a discussion of bifurcations in W 1 ' 1 .) 

The existence of the EAW branch has been inves- 
tigated in many electrostatic Particle-In-Cell Q and 
Vlasov-Poisson 0, [l(| simulations. Moreover, recent 
Vlasov- Yukawa simulations (with Vlasov ions and linear 
adiabatic electrons) led to the prediction of undamped 
electrostatic waves at low frequencies (of the order of 
the proton plasma frequency), similar in nature to the 
EAWs and dubbed ion-bulk waves [HH3]. EAW- type 
fluctuations have been detected in spacecraft data from 
observations in the interplanetary medium (l3j . Also, 
the excitation of the EAWs has been obtained in labo- 
ratory experiments with nonneutral plasmas [Til . [l5j , in 
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which an external driving electric field is applied to the 
plasma column for the time needed to create a popula- 
tion of trapped particles, with the waves surviving after 
the drive is turned off. The flat region in the particle 
velocity distribution generated by trapping at the wave 
phase speed inhibits Landau damping and allows the 
EAWs to survive. The experimental results discussed in 
Refs. [HI, EH confirm the existence of EAWs on the non- 
neutral analog of the so-called thumb curve of Ref. 
as is the case for the numerical results in Refs. [a-TlOj: 
however, these experiments also suggest that wave exci- 
tation can be obtained off of the usual thumb curve of 
Ref. [H, which we will refer to as off-dispersion EAWs. 
The purpose of this paper is to shed some light on the 
nature of these off-dispersion EAWs, by examining the 
sensitivity of the thumb curve to small deviations from 
the Maxwellian used in Ref. Q. 

However, before starting our analysis, we note that in 
a separate research thread, motivated by experiments on 
nonlinear laser plasma interactions [Hj], Afeyan and col- 
laborators carried out extensive Vlasov simulations 
and observed states that they referred to as Kinetic Elec- 
trostatic Electron Nonlinear (KEEN) waves. These large 
amplitude structures also exist off-dispersion; however, 
the KEEN waves observed by these authors are large 
amplitude, while we focus on relatively low amplitude 
EAWs, where any effect of nonlinearity is limited to a 
narrow velocity range of trapped particles. 

Specifically, we investigate both on-dispersion and off- 
dispersion EAWs with a simple linear dispersion analysis. 
As noted above, there is no a priori reason the distribu- 
tion in the vicinity of v$ should be Maxwellian, since the 
structure near v<j, is determined by the history of forma- 
tion of the waves, and since experimentally fine details 
of the distribution function are difficult to measure, it is 
natural to examine the alteration of the thumb curve of 
Ref. 3] caused by small deviations from the Maxwellian. 
Of main importance here is the alteration caused by re- 
placing the trapped particle region of the velocity distri- 
bution by a plateau, but we also consider the alteration 
caused by altering the tail of the distribution function. 

For the plateau distribution function (see Eq. ©) the 
Landau velocity integral is carefully evaluated using a 
high resolution trapezoidal scheme. We look for roots 
of the dispersion function in the high frequency range 
of electron modes, treating the ions as a stationary neu- 
tralizin g backg round charge. Consistent with the exper- 
iments [l5| the analysis shows that undamped EAWs 
exist in a wide range of the (k,oj R ) plane, that is, off the 
thumb curve of 

From the perspective of the linear dispersion analysis, 
the existence of the off-dispersion modes is easy to un- 
derstand. The Landau velocity integral in the dielectric 
function obtains contributions from velocities that are 
well away from the plateau (the non-resonant particle 
contributions) and contributions from near the plateau 
(the resonant particle contributions). The thumb dis- 
persion curve is determined exclusively by contributions 




0.1 0.2 0.3 0.4 0.5 0.6 

k 




1 r i 

I =- 

0.1 0.2 0.3 0.4 0.5 0.6 
k 

FIG. 1: Thumb curve in the k — uu R plane displaying branches 
of undamped LAN waves and EAWs (top). The same thumb 
curve plotted in the k — plane (bottom) . 



from the non- resonant particles, that is, for (k,uj R ) on 
the dispersion curve, the non-resonant contribution alone 
yields a dielectric function that is zero. Off the dispersion 
curve, the non-resonant contribution yields a dielectric 
function that is not zero, so the resonant contribution 
must make up the difference, yielding a total dielectric 
function that is zero. Thus, for the off-dispersion modes, 
electrons in the resonant (or plateau) region make a sig- 
nificant contribution to the mode charge density. In that 
sense, the off-dispersion modes are like beam modes, for 
which a significant part of the charge density resides on 
the beam [17|. As we will see, for the case of the plateau 
(rather than a beam), the resonant particle charge den- 
sity is associated with the two corners of the plateau. 
Thus, we call these waves corner modes. 

As one would expect, the charge density from the cor- 
ners is a very sensitive function of v$ in the plateau re- 
gion. Equivalently, the dielectric function has a spiky 
variation for v$ in the plateau region. Thus, a small 
change in v$ can make the resonant particle charge den- 
sity have whatever value is needed to compensate for the 
non-zero value of the non-resonant dielectric, yielding a 
total dielectric that is zero. On-dispersion EAWs are spe- 
cial only in that the charge density from the two corners 
is equal and opposite adding to zero. This is the case 
if the phase velocity is equidistant from the two corners, 
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that is, at the velocity mid-point of the plateau. From 
this perspective, there is little difference between the 
on-dispersion and off-dispersion EAWs. However, there 
is a significant difference between the EAWs (or corner 
modes) and weakly damped Langmuir waves, where is 
well out on the tail of the velocity distribution and there 
is no significant corner contribution to the mode charge 
density. 

Another way to obtain off-dispersion waves is to alter 
the tail of the distribution function, which may not be 
known precisely. Although the dispersion curve is not 
as sensitive to this kind of deviation, one can ascertain 
systematic shift of roots off of the thumb curve. In par- 
ticular, we can show analytically that a fattening of the 
tail of the distribution shifts roots toward lower /c-values 
and chopping the tail shifts them toward higher k- values. 
Chopping produces perturbed corner charge and this idea 
leads to a derivation of a rule of thumb for assessing shifts 
caused by general plateau type of equilibrium distribu- 
tion functions. 

The paper is organized as follows. In Section II we nu- 
merically analyze the roots of the electrostatic dielectric 
function and discuss the wave dispersion relation for a ve- 
locity distribution flattened in a small velocity interval. 
This is followed by an analysis of the consequences of al- 
tering the tail and the derivation of the rule of thumb. In 
Section III the numerical results of Vlasov-Poisson sim- 
ulations are presented and compared to the analytical 
predictions of Section II for the plateau distributions. 
Summary and Conclusions are given in Section IV. 
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FIG. 2: Velocity dependence of the function f p (top) and its 
first velocity derivative df p /dv (bottom) in the region near 
Vo = 1.50. Dashed line of top plot is the usual Maxwellian, 
f M ; red-vertical lines of the bottom plot mark the edges of 
the interval [V - AV p /2, V + AV p /2]. 



II. WAVE DISPERSION RELATION 



iD I (k,uj R ) + iu) I dD R {k,oj R )/du> R , where 



The propagation of electrostatic waves in a collisionless 
unmagnetized plasmas can be described by a simplified 
1+1+1 (one space, one velocity, and one time dimension) 
Vlasov-Poisson system, which is given in dimensionless 
form as follows: 



df df df dE 



where / = f{x,v,t) is the electron distribution function 
and E = E(x,t) the electric field. In ([TJ, the ions are 
a neutralizing background of constant density no = 1, 
time is scaled by the inverse electron plasma frequency 
ojp , velocities by the electron thermal speed v t h, and 
lengths by the electron Debye length A^- For simplic- 
ity, all the physical quantities will be expressed in these 
characteristic units. 

By linearizing Eqs. ((TJ) and following the Landau pre- 
scription |Il.lT8j for weak wave damping, the time asymp- 
totic solution for the complex frequency of the fluctua- 
tions (uj = uj r + iu) 1 ) can be obtained by looking for the 
roots of the dielectric function D(k,uj) ~ D R (k,uj R ) + 
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D, 
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/oV (2) 



Here, fd(v) := dfo/dv, f indicates an integral over 
all v G K with the singularity handled by taking the 
Cauchy principal value, /o is the equilibrium velocity 
distribution of electrons, and = uj R /k is the wave 
phase speed. The roots of D R give the real part of the 
wave frequency u R , while the imaginary part is given by 
u; I = -D I /(dD R /du; R ). 

Undamped waves can be obtained from Eqs. {3J by 
assuming /o has a velocity plateau of vanishing velocity 
width at v = I),)). This renders D^v^) = and solution 
of D n = yields oj r = u n (k). In Ref. [H, the equation 
D R = was solved by assuming Maxwellian /o with the 
velocity plateau of vanishing width, leading to the so- 
called thumb curve, the dispersion diagram displayed in 
the top plot of Fig. [TJ The upper branch of this k — u n 
diagram represents Langmuir (LAN) waves, while modes 
of the lower branch are usually referred to as EAWs, since 

(EAW) 

for small wavenumbers on the lower branch oj r « 
1.3 lfc, which is reminiscent of acoustic waves. In the 
bottom plot of Fig. [lj the same thumb curve is displayed 
in the k — v^, plane. 
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FIG. 3: Real part of the dielectric function D R as a function 
of vj, for k — 0.3 for a Maxwellian equilibrium velocity distri- 
bution (top) and for the equilibrium distribution f p given in 
Eq. © for Vo = 1 (bottom). 

From the plots in Fig. [T] it is evident that no undamped 
roots of D R exist beyond a critical value of the wavenum- 
ber k* ~ 0.53. The presence of this nose-like structure at 
k ~ k* appears unphysical, since the group velocity of the 
wave at k* seems to diverge. To understand this apparent 
paradox one should bear in mind that the thumb curve 
does not represent a usual dispersion relation: each point 
in the {k,uj R ) plane along the thumb curve corresponds 
to a different particle velocity distribution function. In 
order to get undamped solutions, the location of the in- 
finitesimal plateau in the electron velocity distribution 
must slide along and always fall at v = v^. This corre- 
sponds to changing the shape of the velocity distribution 
at each v^. 

As discussed in Sec. HI the existence of the EAW branch 
has been reproduced in numerical simulations and ob- 
served in experiments on nonneutral plasmas. In the 
experiments, stable oscillations also are observed off- 
dispersion, i.e., off of the usual thumb curve. In order 
to obtain insight for understanding this experimental be- 
havior, we will analyze the roots of D a for an equilib- 
rium particle velocity distribution function that deviates 
from Maxwellian by a small, but not infinitesimal, veloc- 
ity plateau of width AV P located at v = Vo. Specifically, 
we chose the plateau distribution function given by 

where f M = exp (— v 2 /2)/^/2tt is the usual Maxwellian, 
n p is an even integer (here n p — 10), AV P — 0.01, and 
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FIG. 4: A zoom of the bottom plot of Fig.[3]near = Vo = 1. 



N is a normalization constant that deviates slightly from 
unity. It is worth noting that f p is smooth in v with 
derivatives up to order n p that vanish at v — V$. 

The distribution f p and its first derivative df p /dv are 
shown in Fig. [2] in the region near Vo (where Vo = 1.5 
for illustrative purposes); the dashed line in the top plot 
represents the function f M . The red- vertical lines in the 
bottom plot indicate the width of the plateau AV P . It is 
clear from this figure that in the interval [Vq— AV P /2, Vb+ 
AVp/2] the first velocity derivative of f p obtains very 
small values. 

We investigate the possibility of getting undamped (or 
weakly damped) plasma oscillations with v<p in the inter- 
val [Vq — AV p /2, Vo+AVp/2]. In this velocity interval, one 
can numerically calculate the value of D R for a fixed k, 
with the imaginary part of the dielectric function being 
negligible because of the plateau of width AV P at v = Vo 
(see the bottom plot of Fig. [2]). 

To calculate the value of D R , we compute the Cauchy 
principal value of the integral of D R of Eqs. @ on a 
uniform velocity grid with a standard trapezoidal scheme. 
The limits of numerical integration are set to \v m \ = 6, 
with the distribution set to zero outside this interval. To 
smoothly resolve the small plateau of width AV^, = 0.01 
in f p , in this interval we use a large number of grid points 
N v = 12000, so as to have Av = 2\v ma J/N v = 0.001 < 
AVp. 

The two plots of Fig.[3]depict the dependence of D R on 
Vrf,, for a fixed value of the wavenumber k = 0.3. In the 
top plot we show the case where /o is Maxwellian, results 
which are essentially equivalent to the thumb curve of 
Fig. [TJ In the bottom plot of Fig. [3] we show the results 
obtained when the equilibrium distribution is chosen to 
be f p of Eq. ([3]) with a small plateau at v = Vq = 1. The 
curve in the top plot displays two roots, the EAW with 

e w ^ \^ anc i the LAN wave with v£ A ~ 3.86, 
both being undamped since D z is assumed to vanish at 
each t>0. Note, the curve in the bottom plot reveals new 
features: the two roots corresponding to the EAW and 



5 



LAN wave (in agreement with those of the top plot) now 
undergo Landau damping (very strong for the EAW), 
since the location of the velocity plateau in f p is now 
fixed at v = Vo = 1. In addition, D R displays a marked 
spike in the region around v<f, ~ Vo = 1. 

In Fig. U we zoom in for a close-up of the spike region 
around Vq and find several extra roots of D R . Among 
these roots, we focus on the one that falls within the in- 
terval [Vq— AVp/2, Vb+AVp/2] marked by the red-vertical 
lines, since the other roots outside this interval are very 
strongly Landau damped. Now the search for the un- 
damped root of D R and related analysis can be restricted 
to a limited and very small region of the velocity domain, 
and we can increase the velocity resolution by a factor 
of 100 so as to determine more precisely the location of 
the root of D R in the interval [V - AV p /2, V + AV p /2}. 
By doing this, we find that the root of D R is located at 
= ~ 0.9966. We also evaluated the corresponding 
imaginary part of the dielectric function from the second 
of Eqs. (]2j, getting a small value D T ~ 2 x 10 -3 . More- 
over, using ui R = kv^ and w 2 = -kD I /(dD R /dv^)\ v *, 

the ratio Rth '■= I^i/^rI — 10~ 6 , meaning that this is 
an almost undamped solution. It is interesting to point 
out that if one chooses the value of Vq (the location of 
the plateau in f p ) so that for a fixed k it falls exactly 
on the thumb curve of Fig. Q] (i.e., it falls exactly on the 
LAN or EAW branch), the roots of D R are found exactly 
at v$ = Vq and they are completely undamped, since 
(df p /dv) v = v<f>=Vo = =*> Dj = 0. 

In order to establish the domain of parameters for 
which electrostatic waves can exist without being Landau 
damped, by including small plateaus in the equilibrium 
velocity distribution, we calculate the minimum value of 
\D R \ in the velocity interval [V - AVp/2, V a + AVp/2] 
for different values of k and Vq. Then min{|.D R |} = 
corresponds to a root of D R . 

The results for min {|-D B |} are displayed in the k — Vo 
contour plot of Fig. [5l Here, the dark area represents 
the region where min{|D H |} = 0, which corresponds 
to weakly damped solutions, while outside this region 
min{|D Ji |} > 0, which means no solutions exist. The 
red-dashed line represents the thumb curve previously 
shown in the bottom plot of Fig. [TJ For the dark region 
of the contour plot in Fig. O for which the roots of D R 
are in the interval [Vo — AV^/2, Vo+ AVp/2], one can eval- 
uate the ratio Rth to ascertain the importance of Landau 
damping for each solution. The maximum value of this 
ratio in the dark region of the contour plot in Fig. [5] is 
R th := max {\w 7 /tv R \} ~ 6 x 10~ 5 , with \cj t /ui r \ be- 
ing exactly null on the thumb curve (red-dashed line in 
Fig. [5]). Therefore, Fig. [5] shows that almost undamped 
oscillations can be obtained with a small plateau in the 
equilibrium velocity distribution in an unexpectedly wide 
region around the thumb curve. More importantly, un- 
damped solutions can be found well beyond the critical 
wavenumber k* predicted by the thumb curve, suggest- 
ing an avenue for understanding the experimental results 
with nonneutral plasmas discussed in Refs. [l4lll5T|. 
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FIG. 5: Contour plot of min{|D n |} in the k — Vo plane; red- 
dashed line is the thumb curve. 



To complete our analysis, we analyze the perturbed 
distribution function of these undamped oscillations, the 
form of which is given by Sf p = f p /(v - v%) @,[l8j|. We 
first assume k = 0.3 and Vq — 1 (as in Fig. [3]) and v*^ — 
0.9966, which corresponds to a mode that is located off 
the thumb curve (see Fig. [I]). In the top plot of Fig.[6j 6f p 
is plotted as a function of v in the region around Vo = 1. 
The red- vertical lines mark the interval [Vo — AV P /2, Vo + 
AVp/2] and a small spike is seen at the location of the 
pole at v = — 0.9966. In addition, two pronounced 
peaks, visible within a velocity interval I v of width ~ 
0.04 around Vo, correspond to the sharp corners at the 
boundaries of the plateau. For this off-dispersion mode, 
the contributions to Sf p due to the two corners are not 
symmetric, because the wave phase speed is not in the 
center of the velocity plateau, Vq. This means that for 
these off-dispersion modes, when 6f p is integrated over 
I v there will be a net contribution from the corners to 
the charge density. 

The situation is different for modes that fall on the 
thumb curve. The bottom plot of Fig. [6] displays Sf p 
for an on-dispersion mode with k = 0.4 and Vo = 1.561. 
As discussed previously, when the values of k and Vq 
are such that the mode falls on the thumb curve, then 
i>0 exactly equals Vq, the center of the velocity plateau. 
This suppresses the pole in the perturbed distribution 
(which is valid for any distribution with a plateau at 
v = Vo, whose first and second velocity derivatives van- 
ish at v = Vq). Indeed, the Landau pole is not visible 
in the bottom plot of Fig. O and now the contributions 
from the two corners are exactly symmetric (but of op- 
posite sign). Consequently, the peaks corresponding to 
the two corners cancel upon integration over I v , yielding 
a negligible contribution to the charge density. 

In reality, the tail of f p of Eq. ([3]) is composed of a part 
that is algebraic and a part that is Maxwellian. This de- 
viation from Maxwellian only matters for large velocities, 
and it is insignificant for our simulations of Sec. lIIIl where 
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the tail is actually chopped at large velocities. However, 
it does raise the question of how altering the tail affects 
the thumb curve. To address this kind of deviation con- 
sider an equilibrium of the form fo(v) — Nif M + N 2 ftaih 
where TVi and N 2 are yet to be determined normalization 
constants, and assume 



/tail 



ftail 



for v < 
for v > , 



(4) 



where > > 0, and we assume continuity at v*, but 
smoothness is not essential. Normalization requires 



1 = N 1 +N 2 +N 2 A(w.) := N!+N 2 +N 2 / dv (0^1 - /„) 



Then, the thumb curve is given by 

f 



k z = idv 



V<f, 



(5) 



(6) 



as plotted in the bottom panel of Fig. [IJ and its deviated 
form with fo(v) is given by 



k 2 = (N\ + N 2 )M + N 2 T(u*) 



where 



T(v*, v<j>) 



-"tail 



dv 



f 

•I A I 



(7) 



(8) 



Equations (O and (|7|) define a one parameter family of 
deviated thumb curves given by 

k 2 = M{ V<j> ) + N 2 (T{v^ v,) - M{ V4> ) A(u*)) , (9) 

where N 2 determines the fraction of particles in the non- 
Maxwellian part of the tail of the distribution. 

We consider two cases: fat tails and chopped tails. For 
both we suppose v* >> v^, so to good approximation 



T 



Fat Tails : A fat tail is one where <p' tei[l 



(10) 



> /m' which is 



the case if entail is a kappa-distribution function which for 
large v behaves as <^taii ~ c/v a , where c > is a constant. 
Then 



T 



a + 1 v? +1 
and © becomes 

k 2 = M + 



dv 4>- 



1 



tail 



cN 2 



a + 1 



M ■ 



a - 1 v*- 1 ' 
(11) 



(12) 



1 — a 



Therefore, for a > 1, which physically is clearly desired, 
the new contribution is negative and the thumb curve 
moves so as to decrease k 2 . This means fattening the 
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FIG. 6: Velocity dependence of the perturbed distribution 
5f p near v = Vb, for a mode with k = 0.3 and Vo = 1 that 
lies off of the thumb curve (top); for a mode with k — 0.4 and 
Vo = 1.561 lying on the thumb curve (bottom). 



tail with v* >> v<f> shifts the thumb curve upwards in 
the bottom plot of Fig. [TJ 

Chopped Tails: For thin tails the situation is not so clear 
cut. First, to remove the Maxwellian tail we set iVi = 0, 
which amounts to setting /o = iV2 (if (u* — v)f M + H(v — 
f*)</Haii), where H is the Heaviside function that is unity 
for positive argument and zero for negative. For thin 
tails, A would be negative, which would tend to move 
k 2 toward larger values, but the T term is more subtle. 
As an extreme case we will chop the tail at by setting 
</> tai i = 0, giving the derivative of fo a jump discontinuity. 
This jump can be removed, e.g. by interpolation with a 
steep slope, but the results do not change much. With 
the chopped choice 



A«) = - / dvf t 



(13) 



making the last term of © positive, as opposed to the 
case of the fat tail. For T we must evaluate the jump 
using f = N 2 ( - S(v - w.)/ M + H{v, - v)f' M ). With 
some manipulation, we obtain 



k 2 = M - 



No 



2tt 



+ (1 + M) 



dv t 



(14) 
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and it remains to determine which of the 'correction 
terms' within the parentheses dominates. To this end 
we make use of the following inequality : 



Run A 



2e 



< I dve~ 



w* > 0. 



Evidently, k 2 increases if 



1 



2(1 + M) 



> 



and a simple calculation shows this is true if 



v 2 > 



1 



M (1 + M) 



(15) 



(16) 



(17) 



Examination of Fig. 3 of Ref. Q reveals that M > 1 
which means for large i?*, a chopped tail shifts the thumb 
curve so as to increase the values of k 2 . This means 
chopping the tail at v* >> shifts the thumb curve 
downwards in the bottom plot of Fig. [TJ 

One can interpret the chopping of the tail as contribut- 
ing an extreme kind of corner charge to the perturbed 
charge distribution. In closing this section we will use 
this idea to obtain a rule of thumb for explaining fre- 
quency shifts due to plateaus. Consider the following 
plateau distribution with extreme corners: 

fep = N M f M [H(v - v+) + H{v- - «))] 

+N p [H{v-vJ)-H(v-v+))\, (18) 

i.e., with jump discontinuities at v± := Vq ± Ay/2. Nor- 
malization of (|18j) requires 
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FIG. 7: Time evolution of the electric field spectral compo- 
nent Ek(t) (normalized by the maximum driver amplitude 
E DM ) for Run A (top) and Run B (bottom); in both plots 
the red (gray) curve represents the function g(t). 



leading order, and assuming v~ < v$ < v+, to avoid 
Landau damping, produces 



M(V4,) + 



J M P J M 



v+ - v$ 



^0 



(22) 



v + 



1 = N U -N M / dvf M +N p AV. 



Differentiating of (|T8jl gives 

fep = N M f'[H(v-v + ) + H(v. 



>))] 



(19) an expression that displays the two corner charge correc- 
tions, which have opposite signs provided fj~' > N p > 
/W and V- < v < v + . The direction of the shift in k 
depends on which dominates. From Eq. (|22p we obtain 
the following compact rule of thumb: 



(N M f M - N p ) [5(v -v+)-5(v- v-)} (20) 

where the delta function terms represent the corner con- 
tributions. Upon inserting (|20p into setting D R = 0, 
making use of (|19p . and manipulating, we obtain 



M + M\N, 



dv /a-. 



N V AV 



N & 

ly mJ m 



v+ - V,/, 



N M f- ] - N P 



(21) 



where f^> := f M (v±). Expression (|2T|) is valid if f M 
is replaced by any homogeneous equilibrium distribution 
function. Now expanding in AV/Vq << 1, retaining the 



M 



(V - v^) 2 - {AV/2) 2 



(23) 



which makes it very clear how the sign is determined. 
Note that the corners produce a waterbag-like denomina- 
tor, as opposed to beam modes [13], and this contribution 
vanishes for — Vo, in agreement with our discussion 
above. Equation (|23l) can be used in a practical sense: 
even though in this derivation fj^) represents the values 
of the Maxwellian just below and just above the plateau, 
their difference can be viewed as a measure of the total 
corner charge contributions, while AV serves as an ef- 
fective plateau width. Thus, the rule of thumb provides 
a general rule for parameter dependencies of frequency 
shifts, one that should be useful for analyzing experi- 
mental data. 
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III. NUMERICAL SIMULATIONS 

Because the results of Sec. HJ are essentially linear in 
nature, we investigate their resilience by resorting to 
simulations of the nonlinear Vlasov-Poisson system of 
Eqs. ([T]). In particular, we concentrate on modes that 
arise from deviations of the thumb curve, off-dispersion 
modes, caused by the small plateaus. 

Our simulations are performed with an Eulerian 
Vlasov code based on the well-know splitting time ad- 
vance method given in Ref. [2fJ. The phase space do- 
main for the simulations is T> = [0, L] X [— v m ,v x ]. 
Periodic boundary conditions in x are assumed, while 
the electron velocity distribution is set equal to zero 
for \v\ > v max — 6. We investigate disturbances near 
the initial equilibrium of Eq. ([3]), with n p = 10 and 
AV P = 0.01 by applying a drive force. The x-direction is 
discretized with N x — 256 grid points, while w-direction 
with N v = 12000. Our goal is to numerically analyze the 
modes predicted by Fig. [5] 

The plasma is driven by an external electric field that is 
taken to be a sinusoidal traveling wave with phase speed 
i>0 d that exactly matches Vq, the location of the plateau 
of f p . The explicit form of the external field is 

E D (a, t) = g{t)E DM sin (kx -u D t), (24) 

where E DM is the maximum driver amplitude, k = 2w/L 
is the drive wavenumber, which is the maximum wave- 
length that fits in the simulation box, lo d — kv^, D is the 
drive frequency, and g(t) = [1 + (t — r)"/Ar n ]~ 1 is a pro- 
file that determines the ramping up and ramping down of 
the drive. The external electric field is applied directly to 
the electrons by adding E D to E in the Vlasov equation. 
An abrupt turn-on or turn-off of the drive field would 
excite LAN waves and complicate the results. Thus, we 
choose n — 10 so g(t) amounts to a nearly adiabatic turn- 
on and turn-off. The driver amplitude remains near E DM 
for a time interval of order At centered at t = r and it is 
zero for t > t s ~ r + At/2. We will analyze the plasma 
response for many wave periods after the driver has been 
turned off. 

The effect of the driver is to prepare a state (i.e. distri- 
bution function), which is then used as an initial condi- 
tion for the undriven Vlasov-Poisson system of Eqs. ([1}. 
This type of initial condition is an example of those called 
dynamically accessible in [2lTj23| . where they were advo- 
cated and discussed in detail. Ultimately, any perturba- 
tion of a known distribution function within the confines 
of Vlasov-Poisson theory must, in fact, be caused by an 
electric field, since there are no other forces available. 
Thus it is physically very natural to consider such initial 
conditions. Dynamically accessible initial conditions are 
also important because they have a Hamiltonian origin 
and, consequently, preserve phase space constraints. Be- 
cause the perturbed distribution function is obtained by 
evaluating the known state on particle orbits, the per- 
turbed distribution function must have the same level 
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FIG. 8: Resonance peak for Run B; red-vertical line indicates 
the value of the driver phase velocity v<j, . 



set topology as the unperturbed and the areas between 
any level set contours must be preserved. For our simu- 
lations, the dynamically accessible initial conditions used 
amount to evaluating the plateau distribution of Eq. ([3]) 
on the orbits (run backwards) produced by the total elec- 
tric field in the interval < t < i ff. Thus, the initial 
condition for the undriven dynamics that beg ins at t Q s 
is a symplectic rearrangement of f p @, [H, |24| . 

We performed N S i m — 64 simulations for various val- 
ues of the plateau velocity position Vq of f p and the 
wavenumber k, in order to numerically investigate the 
predictions of Fig. [5j For each simulation we chose 
At = 20T, where T = 2ir/kv ( j )D is the wave period 
(that is the plasma is driven for 20 wave periods) and 
t = 80T. Typically the maximum time for the simu- 
lations is t max = 250T, but when needed the system 
evolution is followed up to 1000T. The driver amplitude 
is set for each simulation by adjusting the driver trap- 
ping time T trap := 2ir/yJkE DM to be larger than t max . 
In particular, we set T tra p = 125Ar = 10i max , so that 
trapping does not play a large role in the system evolu- 
tion, i.e. the distribution function changes little during 
the simulation. 

Numerical results for E, show two different kinds of 
electric field response for t > t a g. To see this we plot 
Ek(t), the electric field /c-spectral component, in the 
semi-log plots of Figs. [3 for two different runs denoted 
by A and B with parameters given in Table HI 



Run 


k 


V 4>D 


A 


0.9 


0.3 


B 


0.7 


1.5 


C 


0.4 


1.561 



TABLE I: Relevant parameters for Runs A, B and C 



In the plots of Figs. [71 the electric signals are nor- 
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FIG. 9: Spatial dependence of the electron charge densities 
p (black-dashed line), pi (red-solid line), and p2 (blue-dot- 
dashed line) for Run B (top) and Run C (bottom). 

malized by the corresponding maximum driver amplitude 
E DM and the red (gray) line represents the function g(t). 
As is evident from the plots, in Run A we observe damped 
(more or less exponentially) oscillations after the driver 
has been turned off, consistent with Landau damping, 
while in Run B we observe a stable electric response, 
consistent with plateau suppression. Figure [5] is a semi- 
log plot of the spectral electrostatic energy, obtained by 
the Fourier analysis of E for t > t ff, as a function of u^, 
for the case of Run B, where a stable plasma response is 
recovered at t > t Q ff . This plot reveals that the electric 
field propagates with a phase velocity v^, near the driver 
phase velocity v ( p D (red- vertical line). 

Moreover, for Run B we evaluated the electron charge 
density p = 1 — J fdv at the end of the simulation. This 
undamped mode with the parameters of Table|T]is located 
off of the thumb curve (see Fig. [I]). The black-dashed 
line in the top plot of Fig |9] is the total electron charge 
density, p, obtained by integrating / over the velocity in- 
terval / = [— v max , v max ], the red-solid line is the charge 
density, p\ , obtained by integrating / over the velocity in- 
terval I v = [1.48, 1.52] near Vq, while the blue-dot-dashed 
line is the charge density, p2, obtained by integrating over 
the complement of I v . Clearly, p = pi + P2- From the 
plot it is seen that the contribution to the charge density 
coming from the small interval I v around Vo, which con- 
tains the sharp corners at the boundaries of the plateau, 
is comparable to or even a bit larger than the contribu- 



tion from the rest of the velocity distribution. The same 
calculation has also been performed for a third Run C, 
whose parameters are summarized in table [U This mode 
falls on the thumb curve. Here, as can be seen in the bot- 
tom plot of Fig. [9l the contribution to the charge density 
from the sharp corners within the interval around Vo 
(red-solid line) is significantly smaller than that of the 
rest of the velocity distribution. These runs provide nu- 
merical evidence that supports the predictions of Sec. |TT] 
for the perturbed distribution function summarized in 
Fig. El 

In order to reproduce numerically the predictions dis- 
played in Fig. [5j we analyzed in detail the results of 
N s im — 64 numerical experiments. For each simula- 
tion we evaluated the real part of the frequency lo r and 
the wave damping rate w ; after the external driver was 
turned off. Each simulation is then characterized by cal- 
culating /ui r I := R from the simulation data. To com- 
pare the simulation results with the contour plot of Fig.[5j 
we use the value R th from Sec. [ill as a threshold to di- 
vide our 64 simulations into two classes: Class 1 for which 
R < R t h an d Class 2 for which R > R th . 

These results are summarized in the k — Vo scatter 
plot of Fig. [101 where the simulations of Class 1 are in- 
dicated by black squares, while those of Class 2 by red 
diamonds. In this figure the red-dashed line indicates as 
usual the thumb curve, while the black-solid lines delim- 
its the dark region of the contour plot in Fig. 03 where 
almost undamped roots of the dielectric function have 
been recovered. Figure [TU] clearly shows that the black 
squares fall within the black-solid line, while the red di- 
amonds lie outside this line; thus, results of the analysis 
of Fig. [5] are well corroborated by the simulations. The 
simulations of runs Run A, Run B, and Run C are indi- 
cated by capital letters in Fig. [TUJ Runs A fall outside 
the black-solid lines, Runs B fall inside the dark region 
with black solid-line boundaries, and Runs C are exactly 
on the thumb curve indicated by the red-dashed line. 



IV. SUMMARY AND CONCLUSIONS 

In the present work, roots of the electrostatic dielectric 
function were analyzed when a velocity plateau of small 
but nonvanishing width is present in the equilibrium ve- 
locity distribution of electrons. The numerical solution of 
the Landau integral, performed through a high resolution 
scheme, allowed us to show that quasi-undamped plasma 
oscillations can be obtained off of the thumb curve of 
Ref. 3]. By solving numerically the Landau integral, we 
noted that the presence of the velocity plateau, even of 
very small width, can highly affect the real part of the 
dielectric function, producing marked spikes within the 
velocity interval around the plateau. In a wide region 
of the k — lo r plane, almost undamped roots of the di- 
electric function were obtained. Examination of the per- 
turbed electron distribution function revealed that most 
of the charge density associated with the off-dispersion 
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FIG. 10: k — Vb scatter plot for Class 1 simulations (black 
squares) and Class 2 simulations (red diamonds). The red- 
dashed line is the thumb curve, while the black curves delimit 
the dark region of the contour plot of Fig. [5] 

oscillations comes from the sharp corners at the bound- 
aries of the velocity plateau, and for this reason we called 
these new modes corner modes. A rule of thumb was de- 
rived by assuming infinitely sharp corners, a rule useful 
for gauging how a plateau shifts roots off of the thumb 
curve. 

Next, these analytical predictions were compared with 
the results of Eulerian Vlasov-Poisson simulations with 
high resolution in velocity space. Our simulations were 
initiated by applying a wave-like external electric field 
to drive the plateau distribution of Eq. (|3|) off of equi- 
librium. The external electric field was turned on and 
off adiabatically in such a way to avoid the excitation of 
usual Langmuir waves. Also, the amplitude of the exter- 
nal driver was chosen to be very small so trapping effects 
were minimized. As discussed in Sec. IIII1 the numerical 
results of these nonlinear simulations corroborated the 
linear results of Sec. |TTJ 

Although we spoke of off-dispersion results, it is im- 
portant to note that all of the roots obtained in Sec. [TT] 
are actual linear electrostatic plasma oscillations: while 
LAN waves are approximate time asymptotic states as 
shown by Landau; EAWs, corner modes, stationary in- 
flection point modes, etc. are all exact linear plasma os- 
cillations for specific homogeneous equilibrium distribu- 
tion functions. In essence, what we are really attempting 
is to find a homogeneous equilibrium distribution func- 
tion that best describes a weakly nonlinear theory. Be- 
cause plasmas can exist in states away from thermody- 
namic equilibrium for substantial lengths of time, there 
is no a priori reason to believe the distribution function 
is Maxwellian, particularly if wave-like disturbances are 
excited. Because of the sensitivity of the dispersion rela- 



tion to equilibrium distribution functions, in the tail but 
particularly near v^, we conclude that the thumb curve 
of Ref. [1] is of limited predictive capability. 

Since the original BGK paper Q, it has been under- 
stood that there is nonuniqueness in the construction 
of these nonlinear modes: a given electric field can be 
consistent with a large class of distribution functions. 
It is difficult to pin down the shape of the distribution 
function for the trapped particle population, because the 
formation of this distribution depends on the time his- 
tory. This is true both in experiments and simulations, 
whether the BGK modes evolve out of instability or arise 
by driving the plasma as we have done here. In the case 
of small amplitude disturbances, we make the point that 
this arbitrariness is the same as that in choosing the ap- 
propriate shape of fo(v) near v^. For very small dis- 
turbances, the stationary inflection point modes of @ 
are natural candidates. For larger disturbances, the off- 
dispersion modes of this paper appear to be an attractive 
alternative. 

Thus, the present work bears on the interpretation of 
recent results of nonneutral plasma experiments that pro- 
vide evidence for undamped off-dispersion modes. In the 
experiments of [bH Il5j , any plateau-like structures in the 
electron velocity distribution are created dynamically by 
means of an external driver electric field that traps res- 
onant particles. After the driving process, the plasma 
is strongly inhomogeneous, with the formation of humps 
and depressions in the particle distribution function, pro- 
duced by the nonlinear dynamics triggered by the exter- 
nal field. Thus, one would think that a linear analysis 
might not be relevant; however, by designing a plateau 
and tail for a distribution function that has modes that 
match the experiment, it appears that one can obtain 
a linear theory consistent with some of the experimen- 
tal results, and possibly even infer information about the 
trapped particles. In fact, we note that the rule of thumb 
explains the frequency shifts observed in the experiments 
of (l5j . but continuing with this line of investigation is 
beyond the scope of the present paper, so we conclude 
here. 
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